Visualize fields on a latitude-longitude grid.
Install dependencies
First let's make sure we have all required packages installed.
using Pkg
pkg"add Oceananigans, CairoMakie, Imaginocean, GeoMakie"Construct a test-bed field
Let's plot a field that lives on a latitude-longitude grid.
using OceananigansFirst create a latitude-longitude grid
Nx, Ny, Nz = 180, 120, 2
grid = LatitudeLongitudeGrid(size = (Nx, Ny, Nz),
latitude = (-60, 60),
longitude = (-155, 25),
z = (-1, 0),
topology = (Bounded, Bounded, Bounded))180×120×2 LatitudeLongitudeGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded} on Oceananigans.Architectures.CPU with 3×3×3 halo and with precomputed metrics
├── longitude: Bounded λ ∈ [-155.0, 25.0] regularly spaced with Δλ=1.0
├── latitude: Bounded φ ∈ [-60.0, 60.0] regularly spaced with Δφ=1.0
└── z: Bounded z ∈ [-1.0, 0.0] regularly spaced with Δz=0.5Let's create a field. We choose a field that lives on the faces of the cells but any field should do.
We set the field value to $\sin^2(3λ) \sin(3φ)$ and see how that looks.
field = Field{Face, Face, Center}(grid)
set!(field, (λ, φ, z) -> sind(3λ)^2 * sind(3φ))181×121×2 Field{Oceananigans.Grids.Face, Oceananigans.Grids.Face, Oceananigans.Grids.Center} on Oceananigans.Grids.LatitudeLongitudeGrid on Oceananigans.Architectures.CPU
├── grid: 180×120×2 LatitudeLongitudeGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded} on Oceananigans.Architectures.CPU with 3×3×3 halo and with precomputed metrics
├── boundary conditions: FieldBoundaryConditions
│ └── west: Nothing, east: Nothing, south: Nothing, north: Nothing, bottom: ZeroFlux, top: ZeroFlux, immersed: ZeroFlux
└── data: 187×127×8 OffsetArray(::Array{Float64, 3}, -2:184, -2:124, -2:5) with eltype Float64 with indices -2:184×-2:124×-2:5
└── max=1.0, min=-1.0, mean=1.40974e-182D visualization
We can visualize this field in 2D using a heatmap. Imaginocean.jl adds a method to heatmap! so that it works with Oceananigans.jl fields.
using CairoMakie, Imaginocean
kwargs = (colorrange = (-1, 1), colormap = :balance)
fig = Figure()
ax = Axis(fig[1, 1],
xlabel = "longitude [ᵒ]",
ylabel = "latitude [ᵒ]",
limits = ((-180, 180), (-90, 90)))
heatmap!(ax, field, 1; kwargs...)
figWe can do the same but with a GeoAxis provided by the GeoMakie.jl package that allows us to easily add coastlines or also use various projections.
using GeoMakie
fig = Figure()
ax = GeoAxis(fig[1, 1],
coastlines = true,
lonlims = automatic)
heatmap!(ax, field, 1; kwargs...)
fig